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We present results of a theoretical study of p-Hh and 0-D2 clusters at low temperature (0.5 K < 
T < 3.5 K), based on Path Integral Monte Carlo simulations. Clusters of N < 21 P-H2 molecules 
are nearly entirely superfluid at T < 1 K. For 22 < N < 30, the superfiuid response displays strong 
variations with N, reflecting structural changes that occur on adding or removing even a single 
molecule. Some clusters in this size range display quantum melting, going from solid- to liquid-like 
as T — » 0. Melting is caused by quantum exchanges of molecules. The largest p-Hh cluster for 
which a significant superfluid response is observed comprises N=27 molecules. Evidence of a finite 
superfluid response is presented for 0-D2 clusters of size up to N=14. Magic numbers are observed, 
at which both types of clusters feature pronounced stability. 

PACS numbers: 36.40.-c 67.90. +z, 61.25.Em 



I. INTRODUCTION 

The search for superfluidity (SF) in condensed phases 
of atomic or molecular species other than Helium is the 
subject of a significant investigative effort, both from 
the theoretical and the experimental standpoints. Fluids 
made of light bosonic molecules such as para-hydrogen 
(p-H 2 ) or ori/io-deuterium (0-D2), whose mass is smaller 
than (or equal to) that of helium atoms, have long been 
regarded as potential superfluids jl.]. However, the at- 
tractive interaction between two hydrogen molecules re- 
sults in an equilibrium crystal phase, with a freezing tem- 
perature higher than that at which a superfluid transition 
might take place Q. 

Confinement, and reduction of dimensionality, are 
deemed plausible avenues to stabilize a metastable liq- 
uid phase down to temperatures sufficiently low that a 
superfluid transition may be observed. However, the- 
oretical studies of p-H 2 films adsorbed on various sub- 
strates p| 0, [E] ) as well as in two mathematical dimen- 
sions @, have so far yielded no hint of possible SF. The 
suggestion was made that SF may occur in a strictly two- 
dimensional (2D) p-H 2 fluid embedded in a crystalline 
matrix of Alkali atoms That prediction was actu- 
ally based on numerical simulations of model systems of 
very small size, and there is evidence that the observed 
superfluid signal is merely a finite-size effect [§| . 

On the other hand, sufficiently small clusters can re- 
main "liquid-like" at significantly lower T than the bulk, 
and therefore SF could occur at temperatures that may 
render its observation possible Q. For example, a re- 
cently introduced experimental technique, known as He- 
lium NanoDroplet Isolation spectroscopy (HENDI), al- 
lows one to investigate a single molecular impurity em- 
bedded in clusters comprising from a few, to several thou- 
sand He atoms [l(J. In these experiments, SF of the 
medium surrounding the molecule (i.e., the cluster) can 
be inferred by studying the rotational spectrum of the 
dopant. Specifically, the observed decoupling of the ro- 
tation of the molecule from that of the cluster signals 



the onset of SF in the cluster. Evidence of SF has been 
found at T < 1 K in helium clusters surrounding a lin- 
ear carbonyl sulfide (OCS) impurity, and, with the same 
technique, in clusters of iV=14-16 p-H 2 molecules sur- 
rounding the same impurity; these clusters are enclosed 
in a larger helium droplet [TJ [l2j . The same experiment 
failed to detect superfluid behavior of o-D 2 clusters. 

These observations constitute the starting point for an 
intense theoretical effort, aimed at gaining a deeper un- 
derstanding of SF in clusters of p-H 2 and, possibly, o-D 2 
molecules. Quantum Monte Carlo (QMC) methods ap- 
pear to be the most powerful theoretical tool to attack 
such an interesting problem, and have already been suc- 
cessfully employed to characterize the structure and SF 
of doped 4 He clusters both at finite and zero tempera- 
ture [ill, [HI, HH, [l6| ■ For instance, Path Integral Monte 
Carlo (PIMC) simulations have yielded evidence of SF in 
(p-H 2 )at clusters seeded with an OCS molecule [OCS@(p- 
H 2 )at], for 10 < N < 17, as well as in smaller complexes 
such as OCS@(p-H 2 ) 5 and OCS@(o-D 2 ) 5 [HI, M, E^. 
Reptation Quantum Monte Carlo simulations have pro- 
vided information about the structure of CO@(p-H 2 )at, 
for N ~ 12, and have clarified some features of the ex- 
perimentally measured infrared spectra [2(| Hl| . 

Obviously, of crucial importance to the full characteri- 
zation of the experiments, as well as to the understanding 
of the microscopic origin of SF, remains the study of pris- 
tine (i.e., undoped) clusters. In this sense, hydrogen is 
even richer a playground than helium. For example, p- 
H 2 /o-D 2 clusters provide a unique opportunity to study 
the effect of mass on the superfluid and structural prop- 
erties of such a small systems, because the o-D 2 molec- 
ular mass is twice that of p-H 2 , while the intermolecular 
interaction is very nearly the same. Moreover, because 
the equilibrium phase of hydrogen is crystalline at low T, 
studying clusters of molecular hydrogen offers insight in 
the evolution of the physical properties as the number of 
particles is increased. 

The first theoretical investigation of pure p-H 2 clusters, 
based on PIMC simulations, was published in 1991; it 
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was restricted to three cluster sizes, namely N=13, 18 
and 33. Two clusters were found superfluid and "liquid- 
like", while (^-£[2)33 was observed to be "solid- like" at 

T < 2 K m. 

In a recent article, we carried out a systematic study 
of the physical properties of p-H 2 clusters, as a function 
of both size and temperature (N < 40 and 0.5 K < T 
< 4 K) [23|. Such a study was based on PIMC simula- 
tions making use of a recently introduced worm algorithm 
(WA) (24j. Specifically, we computed superfluid fraction 
ps and density profiles, pointing out the remarkable de- 
pendence on N of ps, for 22 < N < 30 and how, in this 
size range, some clusters display coexistence of insulat- 
ing (solid-like) and superfluid (liquid-like) phases, with 
the latter, induced by zero-point motion and quantum 
exchanges, becoming dominant as T is lowered. We re- 
fer to this peculiar behavior as "quantum melting" . Al- 
though experimental data for pristine clusters are not yet 
available, novel techniques based on Raman spectroscopy 
may soon be able to isolate pure clusters and characterize 
their superfluid properties [251 ]. 

The purpose of the present paper is twofold; on the 
one hand, we provide additional results for P-H2 clusters, 
chiefly the energetics but also for the superfluid fraction 
(ps(T)) and the structure of the clusters as a function 
of their size, also compare results obtained with differ- 
ent intermolecular potentials. We also offer an expanded 
discussion of the material presented in Ref. [23j. More- 
over, we present novel results for pure 0-D2 clusters of 
size up to N=2Q molecules, down to T=0.5 K. We re- 
port superfluid behavior of 0-D2 clusters of as many as 
14 molecules; in these clusters as well, structure and su- 
perfluid response are closely related, and non-trivially 
dependent on N. Specific clusters [i.e., (o-D 2 )i3 and (0- 
02)19] feature greater stability, confirming the existence 
of "magic" numbers N of molecules. 

The reminder of this manuscript is organized as fol- 
lows: in the next section, we describe our microscopic 
model Hamiltonian and offer basic details of our simula- 
tions (the reader is referred to Refs. [H for a full descrip- 
tion of the WA). In Sec IIIII we present and discuss our 
findings. In particular, we provide additional details on 
our results previously published in Ref. and illustrate 
our new numerical data for both P-H2 and 0-D2 pure clus- 
ters. Finally, we outline the main physical conclusions of 
our work. 



II. METHODOLOGY 

We model the system of our interest as a collection 
of N P-H2 (0-D2) molecules, regarded as point particles. 
The quantum mechanical hamiltonian is the following: 

N 

H = -\J2^ + V(R) (1) 



where A=12.031 KA 2 for p-H 2 (half this value for o-D 2 ); 
R = ri, Y2 ■ ■ ■ yn is a collective coordinate referring to all 
N particles in the system, and V(R) the total potential 
energy of the configuration R. 

In most of the published numerical work on clusters of 
molecular hydrogen, the potential energy V is expressed 
as a sum of pair-wise contributions, each represented by 
a spherically symmetric potential; this work is no excep- 
tion, in this regard. For the majority of our calculations, 
we have made use of the Silvera-Goldman (SG) pair po- 
tential [26| . which is the most commonly employed in 
previous investigations. For comparison purposes, how- 
ever, we have also obtained results with a different inter- 
molecular potential, due to Buck [2?], [H, [29[ . 

Our simulations are carried out by means of a canon- 
ical variant of the WA presented in Refs. 24. The only 
difference, with respect to the fully grand canonical im- 
plementation, is that all the so-called Z-sector configu- 
rations (which are utilized to compute thermodynamic 
expectation values) have the same number of particles 
(i.e., there are no fluctuations of N); this proves more 
convenient than a fully grand canonical implementation, 
as we are interested in characterizing properties of clus- 
ters as a function of N. As in any PIMC simulation, a 
high-temperature approximation for the many-body den- 
sity matrix is required, and we have used one that is ac- 
curate up to fourth order (30j in the imaginary time step 
e; we empirically found the value e = 1/640 K _1 to yield 
estimates indistinguishable, within our quoted statistical 
uncertainties, from those extrapolated to the e — > limit 
(i.e., the limit where the method becomes formally ex- 
act). We used in our calculations a cubic box of 50 A 
side, with periodic boundary conditions (PBC) [3l[. I n 
previous calculations (e.g., Ref. |22| ). an artificial confin- 
ing potential was adopted, in order to prevent molecules 
from evaporating (i.e., to keep the clusters together). We 
found this device unnecessary, as all the clusters studied 
here stay together, with no sign of evaporation, without 
such an external potential. 

As mentioned above, in the present work we computed 
energetics, radial density profiles and the superfluid frac- 
tion ps(T) of P-H2 and 0-D2 clusters. All of the estimates 
reported here, as well as the physical effects discussed in 
the next section, are observed to be independent of the 
initial configuration of the clusters. 

In a finite system, and within the PIMC formalism, 
it can be shown [32j that ps(T), defined as the fraction 
of the system that decouples from an externally induced 
rotation, is expressed by: 

4to 2 T , 

Ps(T) — — J— (A ) (2) 
h i c 

Here, (...) stands for thermal average; I c is the classi- 
cal moment of inertia of the cluster, and A the total 
area swept by the many-particle paths, projected onto a 
plane perpendicular to one of the three equivalent rota- 
tion axes. 
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III. RESULTS 

In this section we present our results for pure p-H 2 
and o-D 2 clusters. Some of the results shown for the 
superfluid density of p-H 2 clusters have appeared in Ref. 
|23l |. The results for o-D 2 clusters are presented in the 
last subsection, in which we also point out differences and 
similarities between p-H 2 and o-D 2 clusters. 



A. P-H2 clusters 

1. Energetics 

The total energy per p-H 2 molecule e(N) as a function 
of the cluster size N, at T=l K, computed using the SG 
potential, is shown in Fig[T] Also shown in the figure are 
the estimates from Ref. l22t For comparison purposes, 
we also furnish energetics based on the Buck potential; in 
this case, binding energies are considerably greater than 
those computed with the SG interaction. For example, 
a cluster of 23 molecules is characterized by an energy 
some 2.5 K more negative than that obtained with the 
SG potential. Such a difference can be simply ascribed to 
the deeper attractive well (roughly 2 K at the minimum) 
of the Buck potential. 

The energy values shown in Fig. [T] for the Buck po- 
tential are generally consistent with those reported in a 
recent paper by Guardiola and Navarro [Hj], who com- 
puted them at T=0 by means of Diffusion Monte Carlo 
simulations (3o| . However, for some clusters our energies 
are significantly lower than those of Ref . 134 . For exam- 
ple, for A=23 the value quoted in Ref. [Mis -29.94(2) 
K, whereas our own value at T=l K is —30.513(15) K. 
At present, the reason for such a discrepancy is unknown. 
We note, however, that our data for the Buck potential 
are in agreement with those obtained independently by 
other authors, using T=0 QMC methods [36|,[37j]- 

The energy per p-H 2 molecule is a monotonically de- 
creasing function of the cluster size; the simple formula 



e{N) = A v + A S N~3 + A C N~ 



(3) 



where Ay, As and Ac are the well-known volume, sur- 
face and curvature terms, fails to yield a reasonable fit 
to our data, as well as an accurate extrapolation of the 
bulk energy (~ —90 K), to indicate that the physics of 
the bulk is not yet approached by clusters of size studied 
here. This should be compared with the case of small 4 He 
clusters (of size up to AT=50), for which the same formula 
provides an acceptable fit to the data, even though the 
extrapolated bulk energy differs from the actual value by 
over 0.5 K [H. 

Because in the case of helium, the equilibrium con- 
densed phase is liquid, no qualitative structural change 
occurs as the size of the cluster is increased. For hy- 
drogen, on the other hand, the structure must evolve 
from liquid- to solid-like, presumably in a non-monotonic 




FIG. 1: (color online). Energy per p-Hb molecule versus clus- 
ter size N, at T = 1 K. Estimates have been obtained using 
the SG (circles) and the Buck (boxes) potential. When not 
shown, statistical errors are of the order of, or smaller than 
the symbol size. Solid line is only a guide to the eye. Also 
shown for comparison are results from Ref. (triangles). 




FIG. 2: (color online). Chemical potential of p-lb clusters 
versus cluster size N. When not shown, statistical errors are 
of the order of, or smaller than the symbol size. Solid line is 
only a guide to the eye. Numbers refer to stable ("magic") 
clusters. 



fashion, but rather going through structures of different 
shapes and geometries. It is therefore scarcely surprising 
that a simple formula like ((3]) should not offer equally 
acceptable a fit as for helium clusters. 

Recent T=0 QMC calculations [M [36| have pointed 
out the existence of magic sizes (i.e., iV=13) correspond- 
ing to particularly stable clusters. In order to have a 
clear visual identification of these magic sizes, we plot 
in Fig. [5] the chemical potential fi as a function of the 
cluster size defined as 1351: 



H(N) = E(N - 1) - E(N) 



(4) 



where E(N) is the total energy of a cluster of N 
molecules. 
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FIG. 3: (color online). Temperature dependence of the super- 
fluid fraction ps(T) for a cluster of 18 p-Eb molecules (circles). 
Dotted line is only a guide to the eye. When not shown, sta- 
tistical errors are smaller than the symbol size. Also shown 
for comparison are results from Ref. [22l ] (triangles). 
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FIG. 4: (color online). Superfluid fraction of p-Hb clusters 
versus cluster size N, at T=l K. Estimates have been obtained 
using the SG (circles) and the Buck (boxes) potential. When 
not shown, statistical errors are of the order of, or smaller than 
the symbol size. Solid line is only a guide to the eye. Also 
shown for comparison are results from Ref. [53 (triangles). 



The chemical potential increases monotonically for 
N < 13; for clusters of greater size, its behavior is rather 
irregular. Local maxima are observed for particular val- 
ues of N (i.e., iV=13 and 26). The corresponding clusters 
are characterized by considerable stability, and rather 
compact structure. 

In the following, we show that p-H^ clusters of size 
N < 22 display properties consistent with a liquid-like 
superfluid character; the stability of (p-H2)i3 is, in our 
view, not the sign of the occurrence of a p articular solid- 
like structure, as proposed in Ref. [34j, but rather of 
first shell completion (which does not rule out liquid-like 
behavior) [2(],[25j]. On the other hand, as clarified below, 
the peak of the chemical potential observed at A^=26 
reflects the distinct solid-like behavior of the cluster (p- 
H2)26- 

2. Superfluidity 

Data for the superfluid fraction ps(T) of (p-H2)is are 
shown in Fig [3] As expected, ps(T) is a monotonic de- 
creasing function of T; this cluster is essentially entirely 
superfluid (i.e., p s (T) ~ 1) at T < 1 K. Then, p s (T) 
drops to a value ~ 0.27 in a temperature interval of 1.5 
K and decreases more slowly at higher T. 

A finite system cannot undergo a phase transition in 
a strict sense; therefore, in order to assign a "transition 
temperature" , one must resort to some arbitrary crite- 
rion. Because the notion of superfluid fraction is scarcely 
meaningful when the average number of molecules in the 
superfluid phase is ~ 1, we empirically define our "tran- 
sition temperature" T c as that at which Nps{T c ) ~ 2. 
Hence the estimated transition temperature is T c « 3 K. 

The superfluid fraction computed with the SG poten- 
tial (circles) as a function of the cluster size at T = 1 K. 
is presented in Fig [5] [38j. Clusters of size N < 22 are 



entirely superfluid (or nearly so, at least within the pre- 
cision of our calculations) , and display a temperature de- 
pendence of the superfluid fraction similar to that shown 
in Fig. [3] for N=18. 

An interesting, non-monotonic trend is observed (Fig. 
H]) for the superfluid fraction at low T of clusters of 22 
< N < 30 molecules. The value of ps(T=l K), is close 
to unity for N=22, and drops to a local minimum if a 
single molecule is added (i.e., for N=23). It then rises 
again to approximately 85%, if another molecule is added 
(JV=24), and remains relatively large for 7V=25. The ad- 
dition of another molecule, from iV=25 to ^=26, again 
causes an abrupt drop of ps, to less than 0.1. On adding 
one more molecule, ps becomes again significant (approx- 
imately 25%), but drops sharply once more at N=29; it 
remains small, while still featuring noticeable oscillations, 
for greater TV values. 

Generally speaking, as the number TV of particles in- 
creases, the physics of a cluster ought to approach that of 
the bulk; in its bulk phase at low T, P-H2 is an insulating 
(non-superfluid) crystal, i.e., its character is completely 
different from that of small clusters (i.e., N < 22), which 
remain liquid-like and are entirely superfluid. Data in 
Fig. [4] show how the evolution from liquid- to solid-like 
does not occur continuously, i.e., bulk properties do not 
gradually emerge when N increases. 

In particular, we interpret the peculiar behavior of ps 
observed for 22 < N < 30, as due to alternating super- 
fluid (liquid-like) or insulating (solid-like) character of 
the clusters. It seems plausible that drastic changes of 
the superfluid fraction, occurring upon adding just one 
molecule, ought to be directly connected with structural 
changes. 

In order to illustrate this point, we show in Fig. [5] 
graphical representations of the three-dimensional struc- 
ture of a P-H2 cluster consisting of 25 and 26 molecules 
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FIG. 5: (color online). Three dimensional representation of 
a cluster comprising 25 and 26 p-Hb molecules at a temper- 
ature of 1 K. Even though the information provided by this 
kind of figure is purely qualitative, the structural difference is 
evident: (p-Hzjis is essentially liquid-like (molecules are higly 
delocalized) , while (p-H2)26 is "solid-like" (molecules are more 
localized and clearly distinguishable). 



at T=l K; these figures are produced as explained in 
Ref. [2(| The cluster (p-H2)26 displays remarkably solid- 
like properties. Its structure consists of three rings of 
five molecules, with four other molecules linearly ar- 
ranged along the axes of the rings, while the remaining 
seven molecules form an outer shell. Although their po- 
sition is smeared by zero-point fluctuations, the various 
molecules in the cluster can be clearly identified, indicat- 
ing that they enjoy a fairly high degree of spatial localiza- 
tion; consequently, exchanges among different molecules 
are highly suppressed (though not completely absent, as 
shown in Fig. [6]), and the superfluid response is weak. 

Conversely, the molecules in the cluster (p-H2)25 can- 
not be clearly identified, and the entire system appears 
amorphous. Because of their pronounced derealiza- 
tion, molecules have a strong propensity to be involved 
in quantum exchanges, hence the large superfluid re- 
sponse observed. The physics is reminiscent of that of 
the insulating-superfluid quantum phase transition ob- 
served in lattice models of hard core bosons, where a 
superfluid, non-crystalline phase appears as the system 
is driven away from commensuration (where it is an in- 
sulator) m. 

Figure [5] shows the statistics of exchange cycles involv- 
ing N molecules for the clusters (from bottom to top) 
(p-H 2 )i 3 , (p-H 2 )i 8 , (p-H 2 ) 2 6 and (p-H 2 ) 3 3 at a tempera- 
ture T=l K. Long exchange cycles are known to underlie 
superfluidity. The value of the superfluid fraction of the 
smaller clusters (i.e., N—13 and 18) is close to 1, and 
exchanges occur involving up to N molecules, with sig- 
nificant frequency. On the other hand, in the case of 
(p-H2)26 and (2^2)33 Ps is less then 0.1; remarkably, 
however, exchange cycles involving as many as 23 and 18 
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FIG. 6: (color online). Statistics of permutation cycles in- 
volving N molecules for the clusters (p-H2)i3 (a), (p-EfoJis 
(b), (^H 2 ) 26 (c) (^H 2 ) 33 (d) at T=l K. 



molecules respectively have finite, albeit small, statisti- 
cal weight. This observation is consistent with that made 
in Ref. [39| (based on a PIMC simulation which did not 
explicitly include exchanges). 

The superfluid fraction computed with the Buck po- 
tential (boxes in Fig. 2]) is in excellent agreement with 
that obtained with the SG, for those clusters that are 
liquid- like (i.e., with TV < 22). Quantitative differences 
(i.e., lower values of ps obtained with the Buck potential, 
as shown in Fig. [4J occur for bigger clusters (which gen- 
erally display more solid-like behavior), due to the more 
attractive character of the Buck potential, which results 
in greater molecular localization. However, the observed 
trend is qualitatively identical. 



3. Quantum Melting 

Some clusters in the size range 22 < N < 30, feature 
a fascinating behavior in the temperature range explored 
in this work; we discuss in detail the case of (p-H2)23, for 
which the superfluid fraction takes a local minimum. Fig. 
[7] shows the values of the superfluid fraction and of the 
potential energy per molecule, recorded during a typical 
Monte Carlo run for a cluster of 23 P-H2 molecules at 
T=1.4 K and T=l K. In particular, we show consecutive 
block averages of ps and V, each block consisting of 500 
sweeps [ill ]. 

Despite the large fluctuations affecting the values of 
ps, two different regimes can be easily identified: one 
in which the superfluid fraction is high, with an average 
value close to 1, and the other characterized by low values 
of ps, with an average value close to zero. The poten- 
tial energy, correspondingly, takes on high (low) values 
when ps is large (small). While the behavior of the super- 
fluid fraction is consistent with coexistence of two phases, 
characterized by large and small superfluid response, the 
information given by the potential energy suggests that 
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FIG. 7: (color online). Potential energy per molecule and su- 
perficial fraction observed during a typical Monte Carlo run 
(see text) for a cluster of N=23 p-Hh molecules at T=l K 
(panels a and b) and T=1.4 K (panels c and d). The coexis- 
tence of two phases can be easily recognized since the averages 
of ps and V simultaneously switch between high (liquid-like 
superfluid phase) and low (solid-like insulating phase) values. 
The liquid-like superfluid phase becomes dominant as T is 
lowered. 
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Blocks (500 sweeps) 

FIG. 8: (color online). Potential energy per molecule and su- 
perfluid fraction recorded during a typical Monte Carlo run 
for a cluster of N—18 p-Hh molecules at T=2 K. In this case, 
ps and V simply oscillate around their average values (ap- 
proximately 0.6 for ps) without featuring the clear, simulta- 
neous "jumps" observed in Fig [7] This system is found to 
be liquid-like in the range of temperature considered in this 
work. 



these two phases have also liquid-like and solid-like prop- 
erties. For, the switching of the average value of the po- 
tential energy between two different regimes, separated 
by some ~ 6 K, can be interpreted as due to the sys- 
tem visiting relatively ordered, solid-like, insulating con- 
figurations (characterized by low potential energy), and 
disordered, liquid-like, superfluid ones. The coexistence 
between two distinct (disordered and ordered) phases is 
consistent with the observation made in Ref. [I^] for 
classical clusters. 

On decreasing the temperature from T—IA K (pan- 
els c and d of Fig. [7J to T=l K (panels a and b of the 
same figure), the superfluid (liquid- like) phase (i.e., large 
values of ps and V ) is observed during a greater frac- 
tion of the simulation time, and becomes dominant as T 
approaches K; the cluster therefore "melts" at low T. 
This process is ostensibly induced by zero point motion, 
and the ensuing exchanges of molecules, whose impor- 
tance increases as T is lowered. In this sense, one could 
state that melting is driven by Bose statistics, i.e., it is 
associated to the energy contribution due to quantum 
exchanges. 

On the other hand, when T is sufficiently high, quan- 
tum exchanges are suppressed and the system "freezes" 
in a solid-like structure. Evidence of such quantum melt- 
ing has also been found for (p-H2)27- Other clusters in 
the size range 22 < N < 30 presumably display the same 
physics, albeit at lower T. It is important to stress that 
the behavior illustrated in Fig. [7] is different than that 
of a cluster that is liquid-like, and simply not 100% su- 
perfluid; indeed, in such a case superfluid fraction and 
potential energy merely fluctuate around their average 
values, with no evidence of the switching displayed in 

(c) 




FIG. 9: (color online). Radial density profiles for a P-H2 
cluster with N=18 [ T=0.75 K (boxes) and T=2 K (triangles)] 
and iV=23 [T = 0.75 K (diamonds) and T=2 K (stars)] p-H 2 
molecules. Statistical errors, not shown for clarity, are of the 
order of 5xl0 -4 A -3 or less. 



Fig. [7J This is precisely what we find for clusters of size 
N < 22 , which arc liquid-like at all temperatures, with 
superfluid (normal) component growing at low (high) T 
(see Fig. [SJ. 

Our interpretation is supported by the evolution with 
T of the radial density profiles p{r) of the cluster, defined 
with respect to its center of mass. In Fig. radial den- 
sity profiles for (p-H 2 )23 and (p-H 2 )is at T=0.75 K (dia- 
monds and boxes) as well as at T=2 K (stars and trian- 
gles) are shown. A cluster comprising 23 molecules 
at a temperature T=2 K, displays a two-shell structure 
with one sharp peak at r ~ 2 A, and a second broader 
one at about 5.5 A. As T is lowered to 0.75 K, the first 
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peak becomes significantly broader, and its height de- 
creases. Therefore, molecules in the inner shell are less 
localized, and the cluster is more liquid-like, with greater 
propensity for quantum exchanges both in the first shell, 
as well as between the first and the second shells. For a 
cluster of 18 molecules, on the other hand, density pro- 
files stay the same, as T is lowered from 2 K to 0.75 
K, featuring one particle in the center of mass [signaled 
by the large value of p(r — > 0)] and an outer shell, sep- 
arated by a shallow minimum (indicative of liquid-like 
structure). The value of the superfluid fraction of (p- 
£[2)18 increases from about 0.57 to roughly 1 when T 
goes from 2 K to 0.75 K; however, this change is hardly 
reflected in the cluster structure. Quite differently, the 
increase of ps in a cluster of 23 molecule, observed in 
the same temperature range, reflects a substantial change 
in the structure, due to quantum melting (as shown in 
Fig. [9| . This phenomenon clearly exemplifies the impor- 
tance of quantum effects on the structural properties of 
small systems, which has been also investigated in clus- 
ters of neon, which display moderate quantum character 

EE S3, EE El . 

We conclude by stressing once again that, even though 
the results presented in this section are obtained using 
the SG potential, similar physics is observed if the more 
attractive Buck potential is adopted. The only quanti- 
tative changes are the values of the superfluid fraction, 
which are generally lower than those computed with the 
SG interaction for a given T . However, quantum melting 
is observed in clusters of the same size, even if the poten- 
tial due to Buck is used, i.e., it is broadly independent of 
the details of the interaction. 



B. 0-D2 clusters 

Next, we discuss the properties of o-D 2 clusters; here 
too, we regard o-D 2 molecules as bosons of spin zero [47$ . 
with a mass twice that of p-H 2 molecules. We computed 
energy, superfluid fraction and radial density profiles for 
clusters of size N ranging from 3 to 20 molecules at 0.5 
K < T < 2.0 K. 
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FIG. 10: (color online). Energy per 0-D2 molecules e versus 
cluster size N at a temperature of 0.5 K. Statistical errors are 
smaller than the symbol size. Solid line is only a guide to the 
eye. 



95 
85 
75 
65 



i 55 
45 
35 
25 







19 




13 









4 6 8 10 12 14 16 18 20 
N 

FIG. 11: (color online). Chemical potential of 0-D2 clusters 
versus cluster size N. When not shown, statistical errors are 
of the order of, or smaller than the symbol size. Solid line is 
only a guide to the eye. Numbers refer to stable ("magic") 
clusters. 



confirms the existence of "magic" sizes, corresponding to 
highly stable clusters [i.e., (o-D 2 )i 3 and (o-D 2 )i 9 ] 



1. Energetics 

The total energy per o-D 2 molecule as a function of the 
cluster size, at T=0.5 K is shown in Fig. [TUJ For 7V=3, 
the energy has a value close to —9 K, which smoothly 
decreases with N, reaching a value of roughly —36 K for 
N=12. On adding one more molecule, the energy drops 
to about —39.5 K (iV=13), and decreases more slowly for 
N up to 17. For greater values of N an important drop of 
the energy is still observed at iV=19. As in the previously 
discussed case of p-H 2 clusters, a fit of our data using Eq. 
|3] is difficult, and does not allow an accurate estimate of 
the bulk energy; moreover, the chemical potential, (see 
Fig. [TTj) displaying sharp local maxima for N—13 and 19 



2. Superfluidity 

Fig. [T2l shows ps(T) for a cluster of 11 o-D 2 molecules. 
The monotonically decreasing trend is qualitatively simi- 
lar to that in Fig. [3]for (p-H 2 )ig. However, the superfluid 
fraction of the cluster (o-D 2 )n decreases more rapidly 
with T. In particular, it is approximately 1 at T < 0.5 
K, and suddenly drops to a value ~ 0.39 (corresponding 
to roughly 4 molecules in the superfluid phase) when T 
increases of 0.5 K. At higher temperature, the superfluid 
fraction becomes essentially zero (in the sense explained 
in the previous section) for T ~ 1.4 K, at which the av- 
erage number of molecules in the superfluid phase is less 
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FIG. 12: (color online). Temperature dependence of the su- 
perfluid fraction ps(T) for a clusters of 11 0-D2 molecules 
Dotted line is only a guide to the eye. When not shown, 
statistical errors are smaller than the symbol size. 
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FIG. 13: (color online). Statistics of permutation cycles in- 
volving N molecules for the cluster (o-D2)n at T=0.5 K 
(lower panel) and T=2 K (upper panel). 



than two. It is important to note that, while the cluster 
OH 2 )i 8 is almost 60% superfluid at T=2 K (Fig. [3]), 
(o-D2)n at the same temperature is non-superfiuid, and 
the probability of observing permutation cycles compris- 
ing more than 3 0-D2 molecules is less than 5 x 10~ 4 . 
The statistics of many-particle permutations in the sys- 
tem (o-D2)n at T=0.5 K (lower panel) and 2 K (upper 
panel) shown in Fig. [T3] points out how, at high T, long 
exchange cycles disappear, as clusters behave more clas- 
sically. 

The size dependence of the superfluid fraction at two 
different temperatures (T=0.5 and 2 K) is shown in Fig. 
[H At T=0.5 K, small o-D 2 clusters (N < 9 ) are entirely 
superfluid, while for larger values of N the trend of ps 
is consistent with a noticeable dependence on the cluster 
size, similarly to what observed for P-H2 clusters with 
22< N < 30. A local minimum is observed for 7V=10; 
then the superfliuid fraction grows to a value close to 1 
for N=ll, and remains relatively large at iV=12. 




FIG. 14: (color online). Superfluid fraction of 0-D2 clusters 
versus cluster size -/V, at T=0.5 K (circles) and T = 2 K 
(boxes). When not shown, statistical errors are of the order 
of, or smaller than the symbol size. Solid lines are only guides 
to the eye. 

Large fluctuations observed during lengthy simula- 
tions, render the determination of the value of ps quite 
difficult for clusters greater than 11 0-D2 molecules. Nev- 
ertheless, we can conclude from our data that clusters of 
up to 14 molecules still feature a significant superluid re- 
sponse, which instead is close to zero for larger clusters. 
As expected, at higher temperatures ps becomes progres- 
sively smaller and many-particle permutation cycles are 
almost absent. 

The behavior shown in Fig. indicative of quantum 
melting, has not been observed for o-D 2 clusters, in the 
temperature range explored in this work. We expect it 
to be observable at lower temperatures. 

3. Structure 

Radial density profiles of the clusters (0-02)7, (o-D2)io 
and (o-D2)n at T=0.5 K are shown in Fig. [T5J The clus- 
ter (o-D 2 )7 features a single-shell structure, with a broad 
peak at about 2.5 A, consistent with a high degree of 
delocalization of the molecules. For iV=10, the prob- 
ability of having the center of the system occupied by 
a molecule becomes significant, and the peak is shifted 
to higher distances by roug hly 1 A. The shallow mini- 
mum which appears at r ~ 2 A, suggests that exchanges 
are frequent, and not restricted to molecules in the first 
shell, i.e., they involve also the molecule in the center of 
the cluster. 

Important structural differences characterize (o-D 2 )n; 
as shown by the small value of the radial density at 
the minimum, and by the increased height of the first 
peak in Fig. [15l the molecule in the center of the 
cluster is strongly localized [see also the high value of 
p(r — > 0)], and scarcely participates to multi-particle ex- 
changes. Remarkably, the superfluid fraction of all of 
these systems is large, close to 1 for N=7 and 11. The in- 
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FIG. 15: (color online). Radial density profiles for 0-D2 clus- 
ters with N=7, 10 and 11, at T = 0.5 K. The inset shows the 
radial density profiles of (o-D2)i3 and (p-H2)i3 at the same 
temperature. Statistical errors, omitted for clarity, are of the 
order of 5xl0~ 4 A -3 or less. 



set of Fig. [TBI compares radial density profiles of (o-D2)i3 
and (p-H 2 )i3 at T=0.5 K; the effect of mass difference on 
the structural properties is evident. Indeed, the density 
profile of (p-H.2)i3 reflects the liquid- like nature of such a 
system, whereas in the case of (o-D 2 )i3, structural details 
(i.e., pronounced first shell peak and the small value of 
the radial density at the minimum) are consistent with a 
more orderly (i.e., solid- like) structure. 

It has been already pointed out how local maxima of 
the chemical potential correspond to magic clusters, fea- 
turing greater stability than others (see Fig. [TTj) . Here, 
we discuss in detail the case of (o-D 2 )ig. In Fig. [Pol we 
show the three-dimensional structure of clusters compris- 
ing N=18 and 19 o-D 2 molecules, at T=0.5 K. The clus- 
ter (o-D2)i8 is made of three rings of five molecules, with 
the remaining three molecules linearly arranged along the 
axes of the rings. On adding a single molecule (iV=19), 
which also positions itself along the axes of the rings, the 
structure becomes symmetric with respect to the plane 
of the central ring. Since both (o-D 2 )i8 and (o-D 2 )ig are 
solid- like and insulating (i.e., molecules are clearly distin- 
guishable and the value of ps is less than 0.1), the sharp 
maximum of the chemical potential at iV=19 appears not 
to be related to noticeable structural changes (i.e., from 
more liquid- to solid-like, as observed for p-H 2 when N 
increases from 25 to 26), but rather to the pronounced 
symmetry of this cluster. 



IV. CONCLUSIONS 

We have studied the low temperature physical prop- 
erties of p-H 2 and o-D 2 clusters, of size up to 40 and 
20 molecules respectively, using an accurate numerical 
method. The superfluid fraction of (p-H 2 )jv clusters dis- 
plays a remarkable dependence on the size N for 22 < 




FIG. 16: (color online). Three-dimensional representation of 
clusters comprising 18 and 19 0-D2 molecules, at a temper- 
ature of 0.5 K. Both systems have solid-like properties, as 
molecules enjoy a high degree of localization. However, (o- 
1)2)19 features greater stability, due to its more symmetric 
structure (see text). 



N < 30 pointing out how the bulk phase emerges non- 
monotonically. Coexistence of liquid-like (superfluid) 
and solid-like (insulating) phases has been shown in some 
clusters in this size range; the dominance of the first 
phase at low T indicates that these systems melt as a 
result of the zero point motion and quantum exchanges, 
freezing instead at high T. Superfluid behavior has been 
also found in o-D 2 clusters of size up to N=14. The 
heavier mass of the o-D 2 molecules is responsible for the 
suppression of ps in larger clusters. Some clusters, fea- 
turing specific numbers of particles, are characterized by 
greater stability; the fact that the "magic number" N=13 
is present for both p-H 2 and o-D 2 clusters, suggests that 
its stability is associated with high symmetry, i.e., it is 
related to the potential energy. These structures, how- 
ever, need not necessarily be solid-like, but rather simply 
have filled coordination shells. Given the intense level of 
experimental activity in this field, it is to be hoped that 
some of the predictions made in this work will be soon 
tested experimentally. 
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Note added. Recently we became aware of a 
manuscript by Khairallah et al. [48j], also investigat- 



ing superfluid and structural properties of para- hydrogen 
clusters by path integral Monte Carlo simulations. Al- 
though some of the model interactions utilized therein 
are not the same as in this work, and therefore a direct 
comparison is problematic, their results appear in quali- 
tative agreement with ours. 
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